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Abstract. We study the spatial fluctuations of transient creep deformation of 
materials as a function of time, both by Digital Image Correlation (DIG) measurements 
of paper samples and by numerical simulations of a crystal plasticity or discrete 
dislocation dynamics model. This model has a jamming or yielding phase transition, 
around which power-law or Andrade creep is found. During primary creep, the relative 
strength of the strain rate fluctuations increases with time in both cases - the spatially 
averaged creep rate obeys the Andrade law et ~ while the time dependence of the 

spatial fluctuations of the local creep rates is given by Aet ^ t^°-^. A similar scaling for 
the fluctuations is found in the logarithmic creep regime that is typically observed for 
lower applied stresses. We review briefly some classical theories of Andrade creep from 
the point of view of such spatial fluctuations. We consider these phenomenological, 
time-dependent creep laws in terms of a description based on a non-equilibrium phase 
transition separating evolving and frozen states of the system when the externally 
applied load is varied. Such an interpretation is discussed further by the data collapse 
of the local deformations in the spirit of absorbing state/depinning phase transitions, 
as well as deformation-deformation correlations and the width of the cumulative strain 
distributions. The results are also compared with the order parameter fluctuations 
observed close to the depinning transition of the 2d Linear Interface Model or the 
quenched Edwards- Wilkinson equation. 



PACS numbers: 62.20.Hg, 68.35.Rh, 05.70.Ln, 05.40.-a 



Spatial fluctuations in transient creep deformation 



2 



1. Introduction 

Recent research has highhghted the importance of fluctuations and heterogeneous 
response in the behavior of materials. The central issue is that the "rheological" laws 
that describe what happens to a sample under loading are coarse-grained from the 
microscopic behavior. From this viewpoint, it is then imperative to understand what 
are the origins of various phenomena in this fleld, and to grasp the consequences for e.g. 
materials science. The main emphasis in what follows starts from the phenomenological 
primary creep law, which states that the deformation rate decays in time as a power law, 
as = ei ~ t'^. The particular value 2/3 is coined as Andrade's creep, originating 
from 1910 and its theoretical roots are still under a debate [H [21 13]. 

Crystalline materials deform plastically via the motion of dislocations. Some 
time ago, it was shown that two-dimensional dislocation assemblies exhibit a yielding 
transition [4J. This is a non-equilibrium phase transition separating a state which is 
asymptotically quiescent from an active (yielding) one as the control parameter, external 
stress, is varied. The order parameter - the strain rate - seems to exhibit a second order 
phase transition at the critical yield stress ac- Since it was shown that such simple 
materials science/physics models exhibit critical phenomena the research related to the 
statistical mechanics of crystal plasticity has exploded to various directions including 
coarse-grained models [5], studies of the properties of the transition [6], experimental 
signatures such as crackling noise as acoustic emission [71 18] and in stress-strain curves 
[9J, characterization of the deformation structures [lO], surface patterning [H], and 
three-dimensional studies [T2] . 

However, e.g. Andrade creep is a concept that has been demonstrated in non- 
crystalline materials from rocks to composites [13]. To outline some particular cases, 
the deformation of amorphous metallic glasses has recently become an active fleld of 
its own, and there as well collective phenomena are now being discussed and studied 
[m Unj. in that case, and in the physics of the jamming of granular assemblies 
[T7l [TSl [To] the origins of the intermittent deformation are being studied with concepts 
that are focused upon such materials, such as "shear transformation zones" or the "cage 
effect" [20] . Whatever the microscopic physics, it is an important question of what kind 
of signatures can be actually found in experiments on deformation on the coarse-grained 
scale. In this work, we tackle this issue by studying the spatiotemporal characteristics 
of creep deformation, that is to say the response of material samples to a constant load. 
As noted above, there are indications that collective phenomena might be of decisive 
importance. A short account of the results has been published recently in [2]. 

The questions that are of primary interest are: i) what kind of fluctuations can 
be seen in the creep of experimental samples and idealized dislocation assemblies, by 
computer simulation? ii) do these show universality beyond the models that would be 
appropriate for the material at hand? iii) what kind of correlations and fluctuations 
ensue, for both the "order parameter" (creep rate, for the Andrade creep at least) and 
the integrated order parameter, i.e. the total creep strain? The outcomes for such 
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questions present a challenge for the statistical physics of materials science. There are 
a few theories relevant for the creep deformation of materials. For the particular case 
of dislocation systems, it appears that modified versions of depinning (DP) transitions 
(of elastic manifolds in random media) or absorbing state (AS) phase transitions might 
be in order [5]. The fundamental problem is how to coarse-grain such a theory from 
assemblies of individual topological defects |21]; however the DP/ AS paradigm offers 
suggestions about what to look for in the empirical data whether from simulations of 
a model that is supposed to adhere to such a phase transition picture, or from the 
experiments. 

In what follows, simple paper samples are found to exhibit a primary creep regime 
characterized by the power law decay of the average creep strain rate, et ~ t~° '', 
analogously to the Andrade law valid for many materials ranging from soft metals 
to ordinary paper [H |22]- The spatial variations of the strain rate are studied by 
the digital image correlation method, and are found to be characterized by a power 
law time dependence of the standard deviation of the local strain rates, Aet ~ t~°'^, 
during the Andrade's creep. Thus, the magnitude of the fluctuations of the strain 
rate decays slower with time than the average creep rate, implying that fluctuations 
become more important in relative terms with time. This feature persists in the 
logarithmic creep. Similar behavior is observed in a discrete dislocation dynamics or 
crystal plasticity model, which is known to exhibit a non-equilibrium jamming/yielding 
phase transition at a critical value of the applied external shear stress. These results 
are further compared with simulations of the two-dimensional Linear Interface Model 
(LIM) / quenched Edwards- Wilkinson (qEW) equation close to the depinning transition. 
Such a simple 2d model provides a convenient "benchmark" system with a simple AS/DP 
transition where the relevant phenomena can be studied in a transparent manner. 

The paper is organized as follows: In the next Section we discuss details of 
the numerical simulations of the discrete dislocation dynamics model, as well as the 
experimental setup and measurement methods of the paper deformation experiments. 
In Section [3] we present the results of both the numerical simulations, including also 
simulations of the LIM/qEW case, and experiments. These are then discussed both 
from the point of view of classical theories of Andrade creep in dislocation systems, as 
well as by considering a picture emerging from interpreting the creep deformation as a 
process occurring close to a non-equilibrium second order phase transition. Section H] 
flnishes the paper with conclusions. 

2. Details of simulations and experiments 

2.1. Dislocation model simulations 

The discrete dislocation dynamics (DDD) simulations are performed by means of a two- 
dimensional model of shear deformation with point-like edge dislocations (See [Ij or 
e.g. |23]) . These point dislocations glide under the influence of the local Peach-Koehler 
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Figure 1. A snapshot from the dislocation system. A two-dimensional system of 
discrete edge dislocations subject to an external shear stress a close to the critical 
value (Tc ~ 0.025, showing also the division into boxes of I = 50b. The colors (red and 
blue) indicate the sign of the Burgers vector of the dislocations. Only glide motion 
of dislocations along x direction is considered for simplicity. We study the spatial 
fluctuations of the strain rate by dividing the system into boxes of linear size I, and 
defining the local strain rates over some time interval At, see Equations ([S]) and (|n|). 



forces or shear stresses acting on them. These are superpositions of the apphed external 
shear stress a and the internal stresses as due to the long range anisotropic stress fields 

of all the other dislocations in the system. Here D = yu/27r(l — z/), with fi the shear 
modulus and u the Poisson ratio of the material, b is the magnitude of the Burgers vector 
of the dislocations. Such a system can also be thought to represent a two-dimensional 
cross section of a three-dimensional assembly of straight parallel edge dislocations. Only 
a single slip system is considered (i.e. dislocation glide is allowed in the x-direction only 
within the xy plane of the system), and no dislocation climb is taken into account. The 
equations of motion are taken to be overdamped, such that 



(2) 



where Vn is the velocity of the nth dislocation, Xd refers to the dislocation mobility, s„ is 
the sign of the Burgers vector of the nth dislocation, and a is the external shear stress. 

In the simulations, dimensionless units are used by measuring lengths in units of 
b, times in units of l/{xDb), and stresses in units of D. The stresses are computed by 
imposing periodic boundary conditions in both directions. The equations of motion are 
integrated numerically with an adaptive step size fifth order Runge-Kutta algorithm. 
The system is initially composed of a random arrangement of A^o such dislocations, 
with Burgers vectors b = ifeu^. parallel to the glide direction (with an equal number of 
dislocations with the + and — signs). To mimic dislocation annihilation occurring in 
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real plastically deforming crystals, two dislocations with Burgers vectors of opposite sign 
are removed from the system if their mutual distance is less than 2b. The random initial 
configurations are first let to relax in the absence of externally applied stresses. During 
this initial relaxation, a significant fraction of the dislocations get annihilated, leading 
to a reduction of the dislocation number from the initial Nq = 1600 for L = 2006 to the 
range of = 500 — 600 after the relaxation. During the subsequent dynamics under 
the influence of an external stress a ~ (Tc, only a small amount of further dislocation 
annihilations take place. Notice that the relaxation process described above before 
applying the external stress is essential to obtain the Andrade law: initial conditions 
with randomly positioned dislocations lead to a different, roughly exponential time 
dependence of the strain rate [24j . 

It is known that the DDD model exhibits a jamming or yielding phase transition, 
at a critical value Uc of the shear stress [25]. Below that the activity eventually stops 
as at absorbing state phase transitions typically. Above, in the thermodynamic limit, 
a finite deformation rate exists. By applying a constant external stress, for the setup 
described above, we find ~ 0.025. For an example of a dislocation configuration 
observed close to the jamming threshold, see Figure [H Notice that the dislocations tend 
to form various metastable structures, such as dipoles and walls. 

2.2. Experimental setup 

The experimental setup is shown in Figures [2] and El The lower clamp is connected 
to a set of metal weights. Typical loads were around 40 N which corresponds to a 
stress of 13 MPa. The samples were 100 mm x 30 mm slices cut from standard office 
paper, loaded in the "cross-machine" direction. Relative humidity and temperature 
were kept in constant conditions at 30% and 22 C°, respectively, and the applied load 
was adjusted to achieve different creep times-to-f allure, which generally were between 10 
and 60 minutes. In the experiment, the sample is fixed between steel clamps mounted 
on an aluminum frame and a steel jig is used to ensure the correct alignment of the 
sample between the clamps. The sample is imaged during the experiment with a digital 
camera. A laser distance sensor is used to measure the elongation of the sample together 
with a digital image correlation analysis. 

A set of experiments with a larger magnification was also made. The typical image 
size was about 3 mm x 4 mm and the samples were specially prepared laboratory 
sheets, where 5% of the fibres were treated with colour prior to fabrication or sheet 
making to enhance the contrast and to produce a sharp natural-type pattern. The 
larger magnification makes it possible to reduce the noise level, which is relative to 
digital image dimension, and to obtain the strain rates and their fiuctuations throughout 
the entire experiment. 

The sample was imaged during the experiment with PCO's 1 mega pixel grayscale 
digital camera, SensiCam 370KL0562. The camera has a low thermal noise ratio due to 
the cooled CCD and 12 bit grayscale resolution. The exposure time in measurements 
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Figure 2. The setup for creep experiments. A load is attaclred to tlie lower clamp 
and its movement is controlled by using pneumatic cylinders. A camera was attached 
to the frame, and a laser distance sensor followed the movement of the lower clamp. 
Images were taken during the experiment at rate 0.1 Hz...l Hz. 



was 200 ms. The sampling frequency of the digital imaging is in the range 0.1. ..IHz. 
2.3. Strain field measurement in experiments 

The strain is defined on an evenly spaced grid on an image, which consists of a discrete 
set of points {i,j) with constant spacing Ad. j refers to the strain direction. From the 
DIG one obtains displacements on each point Ay*'-' in a time interval dt. A spatial strain 
rate e^-' in a grid point {i,j) is computed using 

= - — 

where I is the length scale within which the strain is measured. The distribution P(eJ'"') 
ensues, which describes the strain rate evolution during the creep experiment. An 
example of the displacement distribution is shown in Figure HI and the computation of 
the local strain rates is illustrated schematically in Figure El 

The starting point of the DIG algorithm is an elastic image registration algorithm 
[261 1271 128] . The algorithm describes both the image and the deformation using the B- 
spline model. The algorithm finds the deformation function by using the multiresolution 
approach for the minimization, so that the image and the deformation model is refined 
every time the convergence is reached. The criteria for the convergence is the sum of 
squared differences (SSD), i.e. the deformation function is applied to the original image 
and the SSD is computed against the deformed image. The algorithm is described 
in [29]. The knots of the B-splines were defined in an evenly spaced grid, with knot 
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Figure 3. The experimental scenario. From a pair of digital images at a time- 
interval At, the local deformations are extracted at a grid. The strain fluctuations are 
measured via the time-dependent standard deviation, and compared to the mean creep 
rate, here in primary /Andrade creep. In the DDD simulations a similar approach is 
used, by subdividing the system of linear size L into boxes of linear size I in which 
the local strain rates el''' are measured. A digital image of a paper sample on a scale 
of 40 mm^. Superposed is a typical deformation grid for a time difference At of 10 
seconds. The color scale indicates the degree of local creep deformation (blue: small, 
red: large). In the background: the experimental setup. The visible speckle pattern 
has been printed, and designed to have a structure and contrast appropriate for the 
Die method on the scale at hand. Figure reproduced from [5]. 

spacing hx h pixels. The exact algorithm for the deformation computation is described 
in [29]. Knot-spacings, i.e. crates 16x16 to 64x64 pixels were used in the computations. 
One method to estimate the error of the displacement field is by using different crates 
for same image pairs: the difference of deformations was computed using two different 
crates and it was of the order of 0.003 pixels. As a conclusion of these tests in various 
experiments [H |30] we can expect that the algorithm is able to find large continuous 
deformations with an accuracy of about 0.1 mm in the sample-scale imaging. The out- 
of-plane deformations put a limit to the accuracy, since the paper thickness is the order 
of 0.1 pixels in the sample scale experiments. 

Note that the sampling frequency of our imaging is such that we do not expect to 
see the single, microscopic "plastic" or yield events [3lj. That is, our results do not 
reflect directly on the microscopic nature of the dynamics that leads to the fluctuations 
we actually measure. Thus we cannot conclude much on the elementary processes and 
their interactions in space and time, except in the coarse-grained sense. This is actually 
also quite true for the DDD model, though all the dynamics is carried out by the mobile 
dislocations - figuring out the causality of the intermittent dynamics is not easy. One 
should also point out that we measure with the technique localized strain fluctuations, 
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Figure 4. Vector field of absolute displacements at t = 10 s in the creep experiment 
on paper. Vector lengths are scaled so as to just avoid an overlap and thus show the 
relative differences of displacement vectors. Sample time to failure tc = 1549 s. 



which however are not a priori nor a posteriori related to any locahzation in the sense 
of the formation of "strain bands" , or "yield bands" . A separate study should be done 
on the possibility of observing related phenomena in the tertiary creep phase. 

When defining and measuring the fluctuations of the strain-field we study only 
the ?/-component of the strain. This is a simplification since also transverse strains 
are observed. This is seen in Figure H] where we show the vector field of absolute 
displacements taken from the primary creep regime of a sample whose total time to 
failure was = 1549 s. As can be seen, the y-component is clearly the dominating 
component of deformation. 

It is furthermore important to point out that the local strain (rate) fluctuations are 
in no obvious way connected directly to the disordered structure. Figure [3] illustrates 
also partly the variations in the sheet transparency, which is related to the local porosity 
and to the mass per unit area and its fluctuations. These are in turn correlated in a 
non-trivial way with the local elastic and inelastic material response. 

2.4- Measuring acoustic emission and typical creep curve 

The strain from an experiment together with acoustic emission (AE) events is presented 
in the main figure of Figure |5l The inset shows an example signal from an acoustic 
emission timeseries, with the events shown. 
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Figure 5. Typical time vs. strain curve from a creep experiment, with also the AE 
energy shown as a function of time. Red curve shows the total strain of the sample 
and the blue dots indicate acoustic emission events during the experiment, with the 
energy scale indicated on the right. Inset: zoom to the acoustic emission signal, where 
we observe two consecutive events with energies E.^ and E^^i greater than a given 
threshold and a waiting time r^. This figure shows typical creep behaviour of paper, 
with the three phases (primary, secondary, tertiary). It also indicates the relative 
absence of crackling noise, that is the small number of acoustic emission events (AE) 
during the primary (no events) and secondary creep of the sample. 

The AE measurement system consists of a piezoelectric transducer, a rectifying 
amplifier and works with continuous data-acquisition. The time-resolution of 
measurements is 10 fis and the data-acquisition is free of dead time. During the 
experiment we acquire bi-polar acoustic amplitudes by piezocrystal sensor, as a function 
of time. The transducer is attached directly to paper and no coupling agent is used. 
Data acquisition channel has 12-bit resolution and a sampling rate of 312 000/s. The 
transmission time from event origin to sensors is of the order of 5 fis. Acoustic channels 
are first amplified and after that saved to a hard disk. The shape of the AE pulses 
can change and attenuate during transmission, but that should have only a minor effect 
on our analysis. Post-processing (thresholding) the amplitude signal we find a discrete 
sequence of events with corresponding energies Ei and arrival times ti. Figure [5] shows 
that in a typical test, one observes the usual creep phases if the test is allowed to 
continue until final failure. The important point is that actual AE starts only close to 
the final failure of the sample. 

It is interesting to note that there is a large degree of variability (Figure [6]) in the 
typical timescales of the creep process from sample to sample, but that the minimum 
creep rate time tm and the failure time tc are linearly related: a Monkman-Grant -like 
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Figure 6. Time to failure vs. time to minimum of the strain rate from the 
experimental data. Time t„i to the strain rate minimum is proportional to time to 
failure tc- The best-fit proportionality constant indicates tm ^ O.SStc- This is similar 
to the typical Monkman- Grant relation [33]. 



universal behavior exists such that ~ 0.83tc [33|[32]. In other words, there is a degree 
of universahty in the behavior of individual samples since the shape of the creep rate 
curves is such that the life-time is preceded by such a minimum around the transition 
to tertiary creep. The role of the strain fluctuations for the origins of the t^/tc-relation 
should be further investigated. 

3. Results and discussion 

In this section we review the main results, and discuss briefly their interpretation. At 
first we consider arguments for the strain rate fluctuations Aej(t) ~ t~'^. These are in 
spirit related to dislocation based theoretical arguments for the Andrade law presented 
by Cottrell [3l] and earlier by Mott [35]. We conclude that as expected these have 
problems with the fluctuation scaling, as they do not include collective phenomena. 
Then, we link the local creep statistics to scaling theories of integrated order parameter 
fluctuations for absorbing state/depinning phase transitions. This is analyzed both for 
the experimental data and the simulations, including a brief study of the order parameter 
fluctuations exhibited by the LIM/qEW equation close to the depinning transition. 

3.1. Strain rate and fluctuations in the DDD model 

In agreement with previous results [H |36] , we find that the mean strain rate for a ^ ac 
behaves in the relaxation phase as 




(4) 



with 6 ~ 2/3, corresponding to the Andrade/primary creep law. 
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Figure 7. The Andrade's law and the fluctuation scaling for the DDD model for three 
different box sizes to compute the local rates, at a ~ cTc- The early-time cross-over 
in the fluctuations is due to the fixed time-intervals of At — 100 at which they are 
computed. Figure reproduced from [5]. 



We then study the spatial fluctuations of the strain rate by dividing the system 
(see Figure [I]) into boxes of linear size /, and defining the local strain rates over some 
time interval At by 



b 
p 



E 

r„€:hox i,j 



Snb 



Xnit + O.SAt) - Xnit - 0.5At) 
At 



(5) 



The possible choices for the linear box size / are limited by l/^^P < I < L, where p 
is the dislocation density (and thus is the mean dislocation-dislocation distance). 

Therefore, for L = 2006, we consider / = 256, 506 and 1006. There are practical limits 
also for At values that can be used: these should be much smaller than the numerically 
feasible simulation times (of the order of 10"^ — 10^ time units), but also as long as 
possible to be able to meaningfully compare the results with experiments in which two 
consecutive images are separated by a "macroscopic" time. Therefore we chose to use 
At = 100 in the simulations as is visible in Figure [71 

The spatial strain rate fluctuations are quantified by the time dependence of the 
standard deviation of the local strain rates, 



Aei(t) 



E(ei'^(t) - (ei'^(t)))^ 



(6) 



For times t > At, Aet{t) is found to decay like a power law in time, Aet ~ t"''', with 
7 ^ 0.5. Thus, the magnitude of the strain rate fluctuations decays more slowly in 
time than the mean strain rate, implying that the role of the fluctuations becomes 
increasingly important with time. 



3.2. Strain rate and fluctuations in the paper experiments 



Similarly to the results of the DDD simulations, the experimental samples are found to 
exhibit a primary creep regime characterized by a power law decay of the average creep 
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Figure 8. Average creep rate and the standard deviation of the local creep rates from 
the sample scale experiments. Data for a typical set of experiments, I = 60mm and crate 
48 X 48, and SSD = 0.1 pixels. The samples average lifetime is 800-1600 seconds. The 
Andrade-to-logarithmic creep transition appears to take place continuously at 150-200 
seconds. The data is averaged over 16 experiments. The effective Andrade exponent 
is about 0.7, in the range 10 to 50 seconds, and becomes close to one (signalling the 
onset of logarithmic creep) above 150 seconds. Figure reproduced from [2]. 



strain rate, 

(7) 

The spatial variations are found to be characterized by a power law time dependence of 
the standard deviation of the local strain rates, 



where the fluctuation exponent 7 = 0.55 ± 0.1 is found. 

The magnitude of the fluctuations of the strain rate decay slower with time than the 
average creep rate, implying that fluctuations become more important in relative terms 
with time. An average over 16 samples is presented in Figure El and a single example 
is shown in Figure |9l The fluctuation scaling result was independent of the length I 
used. In the sample scale images during creep, the fluctuation scaling was observed 
using / = 14. ..60 mm. Here, the crate is the lower and the sample size the upper bound 
of /. At lower stresses, logarithmic creep is found with no signature of the Andrade 
phase. The DIG analysis on the magnified images shows that the relative strength of 
the fluctuations increase also during the logarithmic creep phase, as evidenced by Figure 
[T0| where / = 0.3 mm was used. 

The average creep rate and the fluctuation amplitude change with time so the 
creep rate probability distribution (PDF) might evolve as well. Typical examples of 
the distributions of (the y-component of) relative strain rates are depicted if Figure [TT| 
showing that such PDF's become narrower but can be otherwise roughly collapsed with 
the expected exponent. 
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Figure 9. An example of the average creep rate and the standard deviation of a single 
sample. The crate is 16x16 pixels, SSD=0.05 pixels and l=3mm. 
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Figure 10. Experimental results from logarithmic creep using a smaller imaging area, 
5 second image intervals, Z=0.3mm, and a smaller load, showing the average creep 
rate and the fluctuations. The former decays faster, and at the end these two become 
comparable, on the observation scale. The data is an average over 9 experiments. 
Figure reproduced from [2]. 



3.3. Local stress-fields and elastic modulus fluctuations of paper samples 

The structure of paper exhibits disorder in many scales [37]. In this work we have 
studied fluctuations of the strain fleld during the creep experiment which corresponds a 
floc-scale and a fiber-scale. The floc-scale corresponds density fluctuations of paper sheet 
at scales order of 1 cm. The flber scale emerges from 0.1 mm thick flbers and flber-flber 
bonds. Density fluctuations lead to fluctuations in the local elastic modulus. When one 
imposes a constant stress during the beginning of the creep experiment, elastic modulus 
fluctuations lead to an inhomogenous stress fleld. The question is can we understand 
fluctuations in the displacement /strain rate fleld based on structural disorder on the 
floe- or flber-scale? 

Andrade's law can be expressed in terms of the global strain e ~ t^^^ or locally 
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Figure 11. The scaled probability distribution functions of strain rates e]'-' from an 
experiment, where the time to failure is 1038 s. The width of the distributions increases 
when time passes. The skewness decreases in this particular sample, but the feature is 
not observed in all experiments. 



el'^i) ~ t^"^^^. The same applies to fluctuations: The fluctuation scaling law Aet{t) ~ 
t~'^ can be studied at different scales, by changing the size of the region where the local 
strain rate is measured, and related to Andrade's law, i.e. to the average strain rate 
over the whole region. We have observed universal behaviour from fiber to sheet scale, 
suggesting that the data needs a model which goes below our observation scale; we do 
not have a similar microscopic description for paper creep as we have for crystallized 
solids based on the concept of dislocation dynamics, i.e. localized topological defects 
that exhibit then collective behavior. 

The inhomogeneous stress field and its relation to the strain rate field is related to 
spatiotemporal aspects of local strain rates during the experiment: when sample creeps 
locally it relaxes the local strain and stress and one would expect that the relative strain 
field is not stationary in time. This is seen qualitatively in accelerated videos of the 
local creep rates: during the Andrade and logarithmic creep, strain fluctuations are not 
localized into fixed regions or zones. 

3.4- Primary creep fluctuations - a simple model 

Here we consider a crude model for such intermittent phenomena. We start with 
the idea of dividing the volume elements/boxes into two categories according to their 
instantaneous activity over the observation time At, i.e. into "active" and "inactive" 
boxes. The active boxes are assumed to be characterized by a local strain rate while 
inactive ones have a strain rate of zero. The probability distribution of the local strain 
rates can then be written as 



P(6h=p5(ei'^-e^) + (l-p)5(eh, 



(9) 
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where p is the probabihty that a randomly chosen box is active. Next we assume a 
power law time dependence for both p{t) and e" (t), 

p{t)^t-", et{t)^t^-' (10) 

such that the mean strain rate obeys Andrade's law, (cj) = (ne") = {n){t1) = Npel ~ 
t~^ . The assumption of two power law -scalings is simply due to the fact that otherwise 
it is difficult to see how to obtain a result similar to the empirical fluctuation results. 
The standard deviation of the local strain rates is then given by 

Aei'^-(t) = {P(t)f''4it)^^-Pit) - t-^'-'^l^\ (11) 

Thus, to get Ae^(t) ~ t"'^ with 7 ^ 0.5, one needs a ~ 1/3. This would also imply 
V ~ Ct'"(^) ~ t^^l"^ ^ i.e. the number of the active boxes and the activity within such boxes 
should have the same scaling with time. Clearly this kind of scaling argument would need 
to be accompanied with a reasoning explaining why the instantaneous activity should be 
proportional to the fraction of active regions/ volume elements. One shortcoming of the 
above crude model is also that dividing different regions into active and inactive ones 
cannot be done unambiguously in systems such as the present ones in which the local 
activity is a (broadly distributed) continuous variable. This makes comparison with 
simulation results with the model difficult, as one needs to threshold the local strain 
rates in order to subdivide the boxes into active and inactive ones. Thus, our attempts 
to verify the above simple model in the simulations of the DDD model or in the paper 
experiments were not conclusive. 

3.5. Classical dislocation-based theories of Andrade creep 

Next we will briefly discuss two classical dislocation-based theories constructed to 
explain the Andrade power law, and demonstrate how such ideas could in principle 
also account for the scaling of the local strain rate fluctuations. Notice however, that 
such classical ideas rely on assumptions that make it difficult to apply them directly to 
the systems of the present study. 

We start by an argument proposed in a study by Cottrell ^34j- There, the strain 
rate is written as a product of the average strain produced by an avalanche, and the 
number of "cages" in the system that initiate such avalanches. By assuming a linear 
work hardening law, these two quantities were argued to be proportional, leading to 
the Andrade creep law. One possibility to account for the observed fluctuation scaling 
would be then to relate these two quantities to e^(t) and p{t) in Equations Q), ffTOl) 
and (fTTj) . respectively, such that scaling of the form of Ael'\t) ~ t^^/^ would follow 
from Equation (fTTI) . Notice, however, that in the theory of Cottrell the avalanches are 
assumed to be initiated by thermal activation, such that an avalanche is triggered when 
a dislocation segment overcomes the pinning force due to a forest dislocation. However, 
neither thermal fluctuations nor forest dislocations are included in the present DDD 
model, and of course in the paper experiment one does not have dislocations as such. 
The underlying implication of the similarity of a dislocation model is that the Andrade 
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law and the fluctuation scaling result in general from localized creep events, which then 
interact via long-range forces. 

Parallel ideas have been earlier presented by Mott |35]: There, a system of 
dislocation sources with the assumption that their activation is due to incoherent 
stress fluctuations by the dislocations that move and thus contribute to the strain was 
considered. One can apply the activation idea directly to individual dislocations, as in 
the present dislocation system (DDD). In this case one would take the dislocations to 
undergo an intermittent burst giving rise to some characteristic strain increment once 
the local stress exceeds some critical value. Together with the assumption of a linear 
hardening law, such an argument was demonstrated to lead to the Andrade law. The 
origin of the linear hardening could then be one of the following: At the critical point 
o" = (Tc of the jamming transition, one could argue that each dislocation moves/relaxes 
(on the average) only once during the experiment, thus leading to a reduction of the 
number of potentially moving dislocations as time goes on. This corresponds in a way 
to an effective "hardening" of the system and gives rise to the time dependence of the 
strain rate. Another more convincing possibility is to assume that the critical values of 
the local stresses to initiate dislocation motion depend on strain. To the lowest order, 
such a dependence would be linear, corresponding to the linear hardening law employed 
by Cottrell and Mott. Again, in this limit, it could be possible to account for the 
fluctuation scaling as in the case of the argument by Cottrell above, by applying ideas 
presented in Equations ffTOj) and f lTTl) . 

As a side note, including higher order terms to the hardening law would then result 
in corrections to the leading Andrade scaling/power-law. E.g. a quadratic law for the 
relation of the stress increment vs. strain increment in a tensile test, Aa ~ OoAe+OiAe"^ , 
with the second RHS term dominating) would imply e(t) ~ t^^^, or an Andrade exponent 
e = 0.8. 

As mentioned above, the classical arguments rely on assumptions (thermal 
activation etc.) that cannot always be justified in the present case. Therefore, we 
shall in the following present ideas based on the general picture of absorbing state phase 
transitions, which we find more appealing for the systems we study here. 

3. 6. Creep fluctuations in an absorbing state phase transition picture 

It has already been suggested that fluctuation phenomena in plastic flow of crystalline 
solids can be described in the theoretical framework of elastic interface depinning [5], 
with an anisotropic long-range elastic kernel (scaling in general as 1/k'm Fourier space). 
Thus, that particular model is very close to the mean-field theory (infinite or high 
dimensional) limit of depinning models. It however fails to describe primary creep as in 
the DDD model, and thus does not serve here. 

Depinning describes the movement of an elastic manifold or interface in the presence 
of a driving force, and quenched or frozen disorder. In the dislocation plasticity 
framework, the idea would be to coarse-grain from discrete dislocations to a continuum 
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field. Then the driving force is the external stress, the disorder is the random stress 
field (from various microscopic origins including large-scale dislocation arrangements), 
and the elastic interactions are coarse-grained from the dislocation ones. A depinning 
transition separates a frozen phase from one with a non-zero order parameter, for 
dislocation assemblies the strain rate, when a critical external stress value is crossed, and 
in the generalized phase diagram the temperature is taken to zero. In this section we 
discuss scaling properties of fluctuations in creep and compare it to what is expected of 
the integrated order parameter at an absorbing state phase transition in full generality. 

The relevant scaling behaviour has been demonstrated using an interfacial (height 
equals integrated activity) description of the contact process (CP) [38j. The contact 
process exhibits an absorbing state phase transition and belongs to the universality 
class of directed percolation. To review its central properties, the CP is usually defined 
in so that each site of the d-dimensional hypercubic lattice is either vacant or occupied 
by a particle. Particles are created at vacant sites at rate A ~ n/2d, where n is the 
number of occupied nearest-neighbors, and are annihilated at unit rate, independent 
of the surrounding configuration. The order parameter is the particle density p and 
the state p = is absorbing. As A is increased above Ac, there is a continuous phase 
transition from the vacuum to an active state. One derives an interfacial model by 
considering the height of the site hi{t) to be the total amount of time that site is 
occupied (or the number of times it has been active). Dickman and Munoz conjecture 
the scaling hypothesis in Equation f|T2|) for the probability density p{h] t) of the height 
at any lattice site at time t [38]. The mean height is {h{t)), and it is expected that 

p{h;t) = jj^^f{h/{hm, (12) 

where / is a scaling function. The mean height then obeys 

m) ~ t'-', (13) 

which implies that the variance of the height Ah"^ scales as 

Ah^ ~ t'-''. (14) 

It should be noted that this kind of scaling is expected in general in models of rough 
interfaces with the simplest case of two independent (spatial and temporal) scaling 
exponents. Thus many usual depinning models mentioned above adhere to it. For 
such models, the early-time roughening exponent P = 1 — 6 (measuring the standard 
deviation of local order parameter variation in time, at a fixed point in space) from 
the above considerations. Moreover as in other cases /3 also describes the one-point 
temporal correlation function in general. 

Fluctuations in the creep deformation of paper and in the DDD can be tested 
against the scaling hypothesis in Equations ( |T2l) and f |T^ . In the experiments on paper 
we measure local relative strains e{x,y;t) at positions {x,y), where local deformations 
are computed as displacements between the initial loaded state and at a time t from the 
initial loaded state. The picture from the initial loaded state is taken 1 second after the 
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Figure 12. Scaling functions of strain fields f{e/{e);t) from experiments on paper 
using I = 14 mm, crate 16x16 and SSD 0.05 pixels. Different lines indicate the time 
at which / is computed, chosen to be at the beginning, in the middle and at the end 
of primary creep. The average distribution over 16 different samples is computed by 
taking all the local displacements at a given time and then computing the scaling 
function from the distribution p{e; t). The scaling functions collapse and are consistent 
with the scaling hypothesis of an absorbing state phase transition (Equation ([T2l) '). The 
distributions have a similar asymmetric form and a large-strain tail which is decaying 
in an exponential fashion. 



load is applied. The scaling hypothesis in Equation ( !T2|) can be applied by interpreting 
the local deformations e{x, y; t) = h{x, y; t) as the local heights of an interface and thus 
the probability density p{e; t) is expected to scale as the one for the local heights in 
Equation ( IT2l) . Notice that here we consider the integrated strain e, and not the strain 
rate e^. 

Equation ( 1131) corresponds to Andrade's law in the creep experiment and Equation 
( IT4|) states that, if the scaling hypothesis holds, the variance Ae^ of the local deformation 
distribution scales as 

Ae^ ~ t^^, (15) 

where 2/3 = 2/3 = 2 — 2^ if the corresponding Andrade's exponent in the Equation f|T3l) 
is taken to be ^ = 2/3. 

In Figure [12] we depict the scaling functions /(e/(e);t) = {e)p{e;t) during the 
primary creep in experiments on paper. The deformation data e{t) is averaged over 
16 samples. Similar result for the scaling of the distribution are obtained for the DDD 
simulations as is shown in Figure [131 The variance Ae^(t) is shown in Figures and 
[T5l During the primary Andrade creep the variance of the fluctuations increases as 
paper experiments and close to t^/^ in the DDD simulations. The 
results indicate that the scaling hypothesis is plausible. The slightly larger exponent 
observed in the case of the paper experiment could be related to the experiments being 
somewhat subcritical. One can also extend the experimental data beyond the range in 
time appropriate for Andrade creep, with at least a rough agreement with the expected 
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Figure 13. Scaling functions f{e/{e),t) from discrete dislocation simulations. 
Distributions of the integrated local strain for three different times for the critical 
stress value CTc ~ 0.025, scaled according to Equation (fT2|) . The box size used is 
I = 256. For larger I (not shown) the scaling function becomes more narrow, possibly 
approaching the shape of the experimental distributions (Figure [T2|) for large The 
both tails of the distributions can be well described by an exponential. 




Figure 14. Variance of the probability distribution of the total deformation p{e; t) 
as a function of time from experiments on paper. The result is the averaged variance 
over 16 experiments. The scaling behaviour corresponds to the Andrade's scaling of 
the experimental data where ej ~ ^-0.7.. .0.8^ This behaviour is in agreement with the 
scaling behaviour shown in the previous figure. 



behaviour for the logarithmic creep phase (an exponent larger than the 0.75 for the 
Andrade regime is found). Notice that the result that the spatial fluctuations of the 
integrated strain exhibit the same scaling as their mean is not in contradiction with the 
observation that the fluctuations of the strain rate scale differently from their mean. 

Finally we analyze the spatial correlation functions of the integrated order 
parameter or the creep strain field. Given that the experimental strain fields are two- 
dimensional, we can study the correlations in two independent directions. The horizontal 
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Figure 15. Variance of the probability distribution of the total deformation as 
a function of simulation time from dislocation simulations, at cr « CTc- The data 
corresponds to that in Figure [T51 The line indicates a scaling derived from an Andrade 
law. A least-squares fit to the data gives a slightly smaller exponent (0.62). The data 
is averaged over 16 samples. 

correlation function is defined as the width 

Wxir, t) = ^ {{e{x + r, - e(x, yc)Y)^, (16) 

where the roughness of the y-directional strain e is averaged over the x-coordinate 
in a constant position Uc- Thus this represents the 1-dimensional correlations of the 
integrated strain for a fixed height Uc which we choose to study here for simplicity, 
given that it is difficult to fix a reference coordinate system on the sample surface. 
Analogously, the vertical correlation function is defined as 

Wy{r, t) = \j {{eiy + r, - e{y, Xc)y)y. (17) 

The scaling behaviour of the width w{r,t) (both and Wy) is expected to be 
governed by the scaling form 

w{r,t)=t'-'f(^^y (18) 

where the above "simple" scaling picture has been assumed to hold, such that P = 1 — 6. 
In Figure [16] we show the squared deformation- deformation correlation functions w^(r, t) 
and Wy{r,t) for three different correlation distances r from a typical experiment. For 
the horizontal direction we observe scaling with an exponent 2 — 26 close to 2/3, 
i.e. consistent with the simple scaling picture and the usual scaling exponent of the 
Andrade's law, 6^2/3. In the vertical direction the statistics is not so good, and 
it is not possible to conclude if the system exhibits anisotropic scaling or not. The 
overall relatively poor statistics also made it impossible to test in detail the possible 
r-dependent saturation as implied by Equation ( |T8|) . 
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Figure 16. An example of the squared deformation-deformation correlation functions 
w{r,t)'^ (Equation (|18|) ') from a single experiment with tc — 6324s, / — 12.3 mm, 
in horizontal (w^(r, i), filled symbols) and vertical {Wy{r,t), open symbols) directions. 
The solid line corresponds to the i^/^ behavior, as derived from the Andrade's law. The 
statistics of the data is insufficient in order to conclude the form of the scaling function 
/(r, t) in Equation (fT8|) . Also attempts to check for the possibility of anisotropic scaling 
failed due to poor statistic in the vertical data. 



3. 7. Spatial order parameter fluctuations at the depinning transition 

For the sake of comparison, we finish by considering the spatial order parameter 
(=velocity) fluctuations exhibited by driven interfaces close to the depinning threshold. 
To this end, we choose to study the d = 2 + 1 Linear Interface Model (LIM), or quenched 
Edwards- Wilkinson (qEW) equation, 

^h^^M. = ^^2h{x, y, t) + r^(x, h) + F, (19) 

where h{x,y,t) is the local interface height at time t, v is the surface tension, r]{x,y,t) 
is a quenched random force term with correlations {r]{x,y,h)ri{x',y',h')) = 2D6{x — 
x')6{y — y')S{h — h'), and F is the external force. We simulate the model in continuous 
time, for linear system sizes up to L = 256. 

The top panel of Figure [T7I shows the spatially averaged interface velocity {v), and 
the spatial fluctuations of the local interface velocities 5v for I = 2 and At = 1 as a 
function of time t for different values of the external force F. Close to a critical value 
F = Fc, both quantities follow asymptotically a power law in time, the average velocity 
obeying (v) ~ t^^, with 6 ^ 0.51. This is in good agreement with earlier results, 
assuming a scaling relation 6 = 1 — P, where /3 ~ 0.48 is the growth exponent [39]. The 
fluctuations follow Av ~ t""^, with 7 ^ 0.25. For F > Fc and F < Fc, both the mean 
and the fluctuations deviate from the power law, approaching a finite constant value for 
F > Fc, and decaying exponentially to zero for F < Fc. 

Considering different length scales / and system sizes L reveals that the growth of 
correlations plays a central role in the scaling of the fluctuations. The characteristic 
power law at F = Fc is observed only after a transient with an / dependent duration. 
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The relative fluctuations Av/ (v) at F = for different / can be collapsed according to 
the scaling form 

Av/{v) = l^'-^'^-'fit/n, (20) 

involving the dynamic exponent z ^ 1.56 [39], and an exponent 6 characterizing the / 
dependence of the fluctuation amplitude. The main figure in the middle panel of Figure 
[T71 shows a data collapse with z = 1.56, 9 — j = 0.26 (~ 6'/2), and 6 = 0.55. This 
collapse suggests that the initial transient is related to the time t* ~ it takes for the 
correlation length ^ to reach the observation scale /. On the other hand, the duration 
of this power law is also limited by the finite system size: The late time data for / = 2 
and At = 1 collapsed according to 

Av/{v) = L<'~^^g{t/L^) (21) 

and shown in the inset of the middle panel of Figure [17] shows that the power law 
describing the relative fluctuations Av/ (v) extends only up to a time t** ~ L^. 

To understand the origin of the fluctuation exponent 7 in this case, it is useful 
to consider the distribution of the local velocities (from which 6v is obtained as the 
standard deviation) as a function of time. The bottom panel of Figure [17] shows that 
the distribution evolves with time from a peaked one towards a power law with a cutoff. 
By normalizing the distributions by the total number of data points (including the 
ones indistinguishable from zero) reveals that the amplitude of the power law decreases 
with time, but otherwise the distribution remains unchanged: Notice in particular how 
the large-velocity cut-off of the distribution is independent of time. Thus, the velocity 
distribution, after the initial transient, can be well described by the superposition 

Piv'^^) = pPm{v''') + (1 - p)5{v'^^), (22) 

where Pm is the distribution of the mobile interface elements, which is found to be 
independent of time after a transient. The full distribution then evolves only via the 
decay of p, describing the fraction of observation boxes with a finite velocity. Thus, in 
the case of the LIM/qEW, an argument similar to the simple box argument. Equation 
( ITTj) . leads to a fluctuation exponent 7 = 9/2 ^ 0.26, in good agreement with the 
numerical estimate 7 ~ 0.25. 

Finally we also consider the (isotropic) height-height correlations 

w{r, t) = sJ{{h{fo + f,t)-h{fo,t)y)r-„ (23) 

analogously to Equations f[T6]) and f[T7j) . After an initial transient, these obey the 
expected scaling form. Equation ( fTSl) . with the scaling function shown in Figure [18] The 
early time growth of the roughness scales with the expected exponent 2/3 = 2—29 = 0.96, 
followed by an r-dependent saturation, controlled by the dynamical exponent z = 1.56. 

The implications of these results for the creep deformation studies are as follows: As 
the fluctuation exponent 7 appears to be larger than 6'/2 in the case of paper experiments 
and the DDD simulations, the box argument. Equation ( ITTi) . would imply that also the 
strain rate distribution of the active elements of the system would evolve in time, to 
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Figure 17. Average and standard deviation of tlie interface velocity in the d = 2 + 1 
LIM/qEW. Top: Average (open symbols) and standard deviation (filled symbols) of 
the local interface velocities for different external forces F in the proximity of the 
critical force Fc ~ 0.25125, for ^ = 2 and At = 1. At the critical point, both exhibit 
a power-law decay in time, with exponents 9 = 0.51 and 7 — 0.25, respectively. The 
fluctuation data has been displaced vertically for clarity. Middle: Scaling of the 
relative interface velocity fluctuations for F = Fc for different I (main figure). The 
inset shows the late time scaling for different system sizes L. In both cases the data 
can be collapsed with the dynamical exponent z = 1.56. Bottom: The distributions 
of the local interface velocities for different times, for F = Fc, 1 = 2 and At = 1. The 
decreasing amplitude in the long time limit implies that an increasing fraction of the 
interface is no longer moving. 
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Figure 18. Scaling of the squared correlation functions w'^{r,t) in the d = 2 + 1 
LIM/qEW in the long time limit, for different correlation distances r. The solid line 
corresponds to 2/3 = 2 — 26* = 0.96, and z = 1.56 controls the r-dependent saturation. 



yield the observed scaling of the fluctuations. Detailed analysis of this kind is however 
difficult as it would require a substantial amount of statistics, and is thus out of scope 
of the present study. Analysis involving growing correlations (in the spirit of Equations 
( l20l) and ( l2Ti) ) could also be complicated by the possibly anisotropic nature of such 
correlations [6]. 

4. Conclusions 

The creep deformation of solids shows, after all, fluctuation features that are unexpected. 
Clearly it is to be expected that in any not completely homogenenous and translationally 
invariant material (sample) the rheology has to include in practice a fluctuating 
component. However, the spatiotemporal aspects that it then encompasses first of 
all are not included in traditional materials science descriptions of empirical data, and 
secondly they may - as indeed is the case here - originate from "collective behavior". 
And, this is the central point of the current investigation. 

The presence of crackling noise in materials is of course an old subject, from 
Barkhausen noise in magnets to acoustic emission in the brittle fracture -related 
deformation of materials and geosystems. In these two particular cases, theory has been 
able to explain to a varying degree from almost complete to qualitative the features 
seen such as probability distributions of noise events, and correlations in the timeseries 
measured. An open problem here, and in many other similar systems where one can 
study collective phenomena in deformation - such that possess a jamming transition for 
instance - either experimentally or numerically or both, is what is the coarse-grained 
theory like, and what is the origin of the eventual criticality and the general picture 
of the fluctuations that ensue. In dislocation plasticity, various attempts of the scaling 
argument -type have been done to understand creep in addition to numerical models 
from DDD to cellular automata. However, no formulation exists that would account for 



Spatial fluctuations in transient creep deformation 



25 



the yielding transition faithfully. Obviously, in non-crystallinc materials the search is 
even further away from a candidate in spite of various rheology-based attempts that have 
been done lately. Thus we conclude that experiments such as published here, including 
the numerical ones, should also offer a valuable guideline for future developments. 

In addition to showing that scaling exists in the fluctuations in primary (Andrade) 
and logarithmic creep, we have also shown the results of the application of various ideas 
from absorbing phase transitions, including depinning transitions for elastic manifolds 
in random media. It appears that the part of creep that would originate from a yielding 
transition agrees with this picture. Our experiments nevertheless suffer from a lack of 
resolution in space and time, and it is to be hoped that a better (two-dimensional) 
test material is found to improve on both. This is true for the creep rate correlation 
functions. Wc note also that the presence of fluctuations in the logarithmic phase is not 
a priori easy to relate to a phase transition picture unlike for the preceding Andrade 
part. The same goes also for the untested behavior in tertiary creep as the sample failure 
is approached. 

Creep is the simplest paradigm of time-dependent failure. What we describe here is 
the development of creep strain mostly when there is not yet actually any damage, i.e. 
in the classical empirical creep language the system is not yet in tertiary creep. This 

overview raises a wide variety of open questions. First, what happens in a ordinary 
tensile test with a certain strain rate, what is the role of the viscoplastic deformation in 
a (apparently) brittle material say, or in a very ductile one if one assumes that similar 
collective phenomena arc to be found there as well? An analogy from depinning is 
a ramp of the driving force so that eventually the threshold is crossed. Second, the 
role of temperature in depinning (and in DDD models) has a specific meaning. It 
infiuences i) the mobihty of the coarse-grained dynamics (scales the timescale in DDD) 
and ii) it leads to "creep" or thermally assisted motion or deformation at the end - but 
beyond the criticality-dependent creep. Third, what is the role of such phenomena in 
materials science, if any? Fourth, one can now go beyond creep to various scenarios of 
time-dependent material failure such as fatigue etc., and repeat the experiments and 
the analysis in such cases. One question is the role of relaxation processes or memory 
effects: how will these influence strain rate fluctuations, and if as for clastic manifolds 
one flnds signatures of the behavior of glassy systems. An example of such could be 
aging, for which many separate scenarios could be envisioned. 

In summary, we have discussed in this article by experiments on a two-dimensional 
material (i.e. paper) and by simulating a crystal plasticity model, that creep deformation 
shows fluctuations that can be analyzed further using the language of absorbing state 
phase transitions. For DDD models the existence of the yielding transition has been 
known, and it has been recently established to be a peculiar, zero-temperature second- 
order phase transition. The observations we have made are of course related to the 
critical scaling thereof: how exactly is still to be understood. The experimental 
signatures from a non-crystalline material are quite similar, which would hint of a greater 
universality than what one might have expected. Future work should include settling 
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some of the open issues listed above, and we would also like to underline the possibilities 
in other systems including model ones such as colloidal particle assemblies, with both 
crystalline order and in the amorphous state ^Ol 1^ . 
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